1 Introduction

This is example illustrates how to apply our new method to calculate the response diversity in a multifarious environmental change context. We will use a dataset of growth rates (our performance trait) of 5 species of phytoplankton factorially grown at 6 different temperatures and 12 phosphate concentrations and their combination in a 1-month experiment. The original dataset is associated with the publication Bestion et al 2018, and the dataset is publicly available on Zenodo. We will apply step by step the procedure described in the main text of the publication and shown with simulated data in the Appendix 1.

In this document we will:

  • Model species response surfaces using GAMs
  • Calculate partial derivatives of one species growth rate to each of the environmental drivers separately
  • Simulate an environmental change scenario that will give us the direction of the environmental change
  • Calculate the response diversity of an hypothetical community assemble by putting together the 5 species of phytoplankton
  • Calculate Response Capacity of two hypothetical communities of 3 species each, assemble by randomly selecting 3 species from the 5 species of phytoplankton

Some of the figures generated in this document are also used in the main text of the publication. To reproduce all the figures shown in the main text, please refer to the R script called “Figure_empirical.R” available in the GitHub repository

1.0.1 Loding data and plotting species responses

Species responses to the environmental drivers. (a). Species responses to temperature (°C). (b). Species responses to light intensity

Figure 1.1: Species responses to the environmental drivers. (a). Species responses to temperature (°C). (b). Species responses to light intensity

1.0.2 Fittig GAMs on empirical data

SpeciesName E1 E2 predicted
Ankistrodesmus 15 0.01 0.6709741
Ankistrodesmus 15 0.51 0.6505641
Ankistrodesmus 15 1.01 0.6501863
Ankistrodesmus 15 1.51 0.6650670
Ankistrodesmus 15 2.01 0.6845927
Ankistrodesmus 15 2.51 0.7001079

1.0.3 Plotting surface for a single species

Figure 1.2: Response surface fitted with GAM. High non-linearity.

2 Partial derivatives for a single species

2.1 E1 - Temperature

First, we calculate the partial derivative with respect to the first environmental variable - temperature - keeping phosphate concentration constant at 15 μmol L−1.

Response surface of Ankistrodesmus The solid line shows at which level of temperature partial derivativ is going to be calculated.

Figure 2.1: Response surface of Ankistrodesmus The solid line shows at which level of temperature partial derivativ is going to be calculated.

Partial effect of temperature on the growth rate of Ankistrodesmus.

Partial effect of temperature on the growth rate of Ankistrodesmus when phosphate concentration is held constant at 15μmol L-1.

Figure 2.2: Partial effect of temperature on the growth rate of Ankistrodesmus when phosphate concentration is held constant at 15μmol L-1.

Partial derivative with respect to temperature on the growth rate of Ankistrodesmus.

Partial derivative with respect to temperature on the growth rate of Ankistrodesmus when phosphate concentration is held constant at 15μmol L-1

Figure 2.3: Partial derivative with respect to temperature on the growth rate of Ankistrodesmus when phosphate concentration is held constant at 15μmol L-1

2.2 E2 - Phosphate concentration

Second, we calculate the partial derivative with respect to Phosphate concentration keeping temperature constant at 27.5 °C.

Response surface of Ankistrodesmus.
Response surface of Ankistrodesmus. The solid line shows at which level of phosphate the partial derivative is going to be calculated.

Figure 2.4: Response surface of Ankistrodesmus. The solid line shows at which level of phosphate the partial derivative is going to be calculated.

Partial effect of phosphate concentration on the growth rate of Ankistrodesmus

Partial effect of phosphate concentration on the growth rate of Ankistrodesmus when temperature is held constant at 27.5°C.

Figure 2.5: Partial effect of phosphate concentration on the growth rate of Ankistrodesmus when temperature is held constant at 27.5°C.

Partial derivative of Ankistrodesmus growth rate with respect to phosphate concentration

Partial derivative with respect to phosphate concentration when temperature is constant at 27.5°C.

Figure 2.6: Partial derivative with respect to phosphate concentration when temperature is constant at 27.5°C.

2.2.1 Plot surface and partial derivatives

Plot the two partial derivatives and relative effects

Summary plot of Ankistrodesmus (a) response surface of Ankistrodesmus (b) Partial effect of temperature on the growth rate of Ankistrodesmus when phosphate concentration is held constant at 15 (μmol L−1). (c) Partial derivative with respect to temperature when phosphate concentration is held constant at 15 (μmol L−1). (d) Partial effect of phosphate concentration on the growth rate of Ankistrodesmus when temperature is held constant at 27.5°C (e) Partial derivative with respect to phosphate concentration when temperature is held constant at 27.5°C.

Figure 2.7: Summary plot of Ankistrodesmus (a) response surface of Ankistrodesmus (b) Partial effect of temperature on the growth rate of Ankistrodesmus when phosphate concentration is held constant at 15 (μmol L−1). (c) Partial derivative with respect to temperature when phosphate concentration is held constant at 15 (μmol L−1). (d) Partial effect of phosphate concentration on the growth rate of Ankistrodesmus when temperature is held constant at 27.5°C (e) Partial derivative with respect to phosphate concentration when temperature is held constant at 27.5°C.

3 Directional deriviatives

To calculate the directional derivatives for all spp used in the experiment, we first create a time series with phosphate concentration and temperature changing randomly over time, which gives us the direction of the environmental change, and then we calculate partial derivatives.

Time series of temperature and phosphate concentration and changing over time.

Time series of (a) temperature and (b) phosphate concentration changing over time.

Figure 3.1: Time series of (a) temperature and (b) phosphate concentration changing over time.

Table with calculated partial derivatives for each sp at different times (only first 6 rows shown).

sp time E1_ref E2_ref pd_E1 pd_E2
Ankistrodesmus 1 33.5 37.51 -0.1354749 0.0080680
Ankistrodesmus 2 26.0 37.01 0.0428851 0.0330899
Ankistrodesmus 3 40.0 29.01 -0.1553823 -0.0000048
Ankistrodesmus 4 24.0 27.51 0.0591728 0.0305064
Ankistrodesmus 5 25.0 26.51 0.0528198 0.0305497
Ankistrodesmus 6 30.5 34.51 -0.0776470 0.0212365

3.1 Calculating response diversity for a community (all 5 spp used in the experiment)

First, we need to calculate the directional derivatives in the trajectory of the env change.

sp time E1_ref E2_ref pd_E1 pd_E2 nxt_value_E1 nxt_value_E2 del_E1 del_E2 unit_vec_mag uv_E1 uv_E2 dir_deriv
Ankistrodesmus 1 33.5 37.51 -0.1354749 0.0080680 26.0 37.01 -7.5 -0.5 7.516648 -0.9977852 -0.0665190 0.1346382
Ankistrodesmus 2 26.0 37.01 0.0428851 0.0330899 40.0 29.01 14.0 -8.0 16.124516 0.8682431 -0.4961389 0.0208176
Ankistrodesmus 3 40.0 29.01 -0.1553823 -0.0000048 24.0 27.51 -16.0 -1.5 16.070159 -0.9956342 -0.0933407 0.1547044
Ankistrodesmus 4 24.0 27.51 0.0591728 0.0305064 25.0 26.51 1.0 -1.0 1.414214 0.7071068 -0.7071068 0.0202702
Ankistrodesmus 5 25.0 26.51 0.0528198 0.0305497 30.5 34.51 5.5 8.0 9.708244 0.5665288 0.8240419 0.0550981
Ankistrodesmus 6 30.5 34.51 -0.0776470 0.0212365 23.5 11.51 -7.0 -23.0 24.041631 -0.2911616 -0.9566739 0.0022915

Then we can calculate response diversity for an hypothetical community containing all the species tested i this experiment.

time E1_ref E2_ref Ankistrodesmus Chlamydomonas Monoraphidium Raphidocelis Scenedesmus rdiv sign Med Med_sing
1 33.5 37.51 0.1346382 0.1872009 0.0864774 0.1041407 0.0282515 1.059515 0.0000000 1.039969 0.3221012
2 26.0 37.01 0.0208176 -0.0195809 0.0446617 -0.0196109 -0.0242417 1.028641 0.7036438 1.039969 0.3221012
3 40.0 29.01 0.1547044 0.2382443 0.4808502 0.1418676 0.4178197 1.154187 0.0000000 1.039969 0.3221012
4 24.0 27.51 0.0202702 -0.0377794 -0.0103271 0.0091463 -0.0520521 1.030865 0.5605527 1.039969 0.3221012
5 25.0 26.51 0.0550981 0.0396279 0.0372339 0.0011792 -0.0102725 1.027214 0.3142855 1.039969 0.3221012
6 30.5 34.51 0.0022915 -0.0042571 -0.0460690 0.0084282 -0.0147952 1.020270 0.3093090 1.039969 0.3221012

Plot response diversity over time

Directional derivatives and response diversity with known trajectory of env change. a. Species directional derivatives over time. b. Response diversity measured as similarity-based diversity metric. c. Response diversity measured as divergence (sign sensitive).

Figure 3.2: Directional derivatives and response diversity with known trajectory of env change. a. Species directional derivatives over time. b. Response diversity measured as similarity-based diversity metric. c. Response diversity measured as divergence (sign sensitive).

4 Response capacity

4.1 Community 1 - Randomly selecting 3 species

Response surfaces of the species composing community 1

Figure 4.1: Response surfaces of the species composing community 1

4.2 Community 2 - Randomly selecting 3 species

Response surfaces of the species composing community 2

Figure 4.2: Response surfaces of the species composing community 2

Response Capacity of (a) community 1 and (b) community 2.

Figure 4.3: Response Capacity of (a) community 1 and (b) community 2.

5 Assessing uncertainty in GAM predictions

To assess the uncertainty associated with GAM predictions and determine the level of replication needed for robustness, one can use posterior sampling. Instead of resampling data (e.g. bootstrapping), one can sample from the posterior distribution of the model. One can either sample new fitted values or one can use the new data (predictions). There are functions in the package gratia to do this: fitted_samples(), predicted_samples(), and posterior_samples(), depending on what one is trying to compute or show the uncertainty of. Here we want to show the uncertainty in the estimated GAM; for that fitted_samples() is sufficient. We will then plot the predicted growth rate with the uncertainty, and compare those with the real / measured data, which help understand the variability in predictions and provide a qualitative assessment of the robustness of the GAMs. More on posterior sampling can be found here: https://gavinsimpson.github.io/gratia/articles/posterior-simulation.html

We first show how it work for a specific subset of the model (i.e. model effect of temperature on one species for a specific phosphate concentration). Then, we assess the uncertainty in the estimated GAMs for all species, all phosphate concentrations, and all temperature.

We use the function “fitted_values()”. Posterior fitted values are draws from the posterior distribution of the mean or expected value of the response. These expectations are what is returned when you use predict() on an estimated GAM, except fitted_samples() includes the uncertainty in the estimated model coefficients, whereas predict() just uses the estimated coefficients.

One species.
Predicted growth rate with the uncertainty vs Measurted values. Black dots are the measured values, and the blue line is the smooth on real data. Adding the posterior fitted samples to the plot of the data, superimposing the Bayesian credible interval on the fitted values we see the posterior draws are largely contained the credible interval.

(#fig:Uncert_sp)Predicted growth rate with the uncertainty vs Measurted values. Black dots are the measured values, and the blue line is the smooth on real data. Adding the posterior fitted samples to the plot of the data, superimposing the Bayesian credible interval on the fitted values we see the posterior draws are largely contained the credible interval.

All species.

Predicted growth rate with the uncertainty vs Measured values. Respose curves along the temperature gradient.

Figure 5.1: Predicted growth rate with the uncertainty vs Measured values. Respose curves along the temperature gradient.

Predicted growth rate with the uncertainty vs Measured values. Respose curves along the phosphate gradient.

Figure 5.2: Predicted growth rate with the uncertainty vs Measured values. Respose curves along the phosphate gradient.

We can see that the posterior draws are largely contained the credible interval.